Preparing the data

Loading the packages and data directories

options(width=108)
out=F
# The list of subjects, the order of conditions, and the thresholds are derived from Subjects.xlsx
library(xlsx)
library(ggplot2)
library(plyr)
library(reshape)
library(matrixStats)
library(gridExtra)
library(lme4)
library(lmerTest)
library(BayesFactor)
#library(splines)
db <- '/home/egor/Dropbox/' # on Linux
db <- '/Users/Egor/Dropbox/' # Windows
# db <- '~/Dropbox/' # on Mac
source(paste(db, 'Prog/R/myFunctions/pvalfn.R', sep=''))

Plot variables

# theme for plotting:
alpha <- .7
w <- .56 # proportion width in group plots
themefy <- function(p) {
    p <- p + theme_bw() + 
         theme(panel.grid.minor.x=element_blank(), panel.grid.minor.y=element_blank(),
            axis.text=element_text(size=8), axis.title=element_text(size=9),
            legend.text=element_text(size=8), legend.title=element_text(size=9),
            legend.key = element_blank(), legend.margin=margin(t=-.04, unit='in'),
            legend.background = element_rect(fill='transparent'),
            plot.title=element_text(face='bold'))
}
cc <- c('#e31a1c','#fdbf6f','#a6cee3','#1f78b4','#b2df8a','#33a02c','#fb9a99')
xLab <- 'Target Peak Time (s)'
yLab <- 'Log Contrast Threshold'
# colLab <- expression(paste('\nTarget\nVelocity (', degree, '/s)', sep=''))
colLab <- expression(paste('Target Eccentricity (', degree, ')', sep=''))
# colLabType <- 'Mask Type'
yLim <- c(-1.5,-0.5)
dodge <- position_dodge(width=0.1)
sumFn <- function(ss, subjStr='subj', xStr='targTpeak', grpStr='targEcc'){
    sumSubj <- ddply(ss, c(subjStr, xStr, grpStr), summarise,
                     mnS=mean(threshMean), se=sd(threshMean)/sqrt(length(threshMean))) #, .drop=F)
    sumSubjMn <- ddply(ss, c(subjStr), summarise, mnStot=mean(threshMean)) # total mean across conditions per subj
    sumSubj <- merge(sumSubj, sumSubjMn)
    sumSubj$normS <- - sumSubj$mnS / sumSubj$mnStot # normalized subject mean
    sumSubj$seNorm <- NA
    # sumSubj$mnS[is.na(sumSubj$mnS)] <- 0
    sumGrp <- ddply(sumSubj, c(xStr, grpStr), summarise,
                  mn=mean(mnS), se=sd(mnS)/sqrt(length(mnS)),
                  norm=mean(normS), seNorm=sd(normS)/sqrt(length(normS)))
    sumGrp$subj <- 'average'
    sumSubj <- rename(sumSubj, c(mnS='mn',normS='norm'))
    sumComb <- rbind(sumGrp, subset(sumSubj, select=-mnStot))
    sumComb$se[is.na(sumComb$se)] <- 0
    sumComb
}
plotAve <- function(pss, subjStr='subj', xStr='targTpeak', grpStr='targEcc', 
                    xlab=xLab, ylab=yLab, collab=colLab, yStr='mn', seStr='se'){
    pss$yMin <- pss[,yStr] - pss[,seStr]
    pss$yMax <- pss[,yStr] + pss[,seStr]
    pss[,grpStr] <- factor(pss[,grpStr])
    p <- ggplot(pss, aes_string(x=xStr, y=yStr, colour=grpStr, group=grpStr,
                            ymin='yMin', ymax='yMax')) +
        geom_point(position=dodge, size=1, alpha=alpha) + geom_line(position=dodge, alpha=alpha) +
        scale_x_continuous(breaks=c(.5,1,1.5), labels=c('0.5','1','1.5')) +
        geom_linerange(position=dodge, show.legend=F, alpha=alpha) +
        labs(x=xlab, y=ylab, colour=collab) + #ylim(0,1) + 
        guides(colour=guide_legend(keyheight=.3, default.unit='inch'))
    p <- themefy(p)
}
plotIndiv <- function(pss, subjStr='subj', xStr='targTpeak', grpStr='targEcc', 
                    xlab=xLab, ylab=yLab, collab=colLab, yStr='mn', seStr='se'){
    pss$yMin <- pss[,yStr] - pss[,seStr]
    pss$yMax <- pss[,yStr] + pss[,seStr]
    pss[,grpStr] <- factor(pss[,grpStr])
    p <- ggplot(pss, aes_string(x=xStr, y=yStr, colour=grpStr, group=grpStr,
                            ymin='yMin', ymax='yMax')) + 
        facet_wrap( ~ subj, ncol=4) +
        geom_point(position=dodge, size=1, alpha=alpha) + geom_line(position=dodge, alpha=alpha) +
        geom_linerange(position=dodge, show.legend=F, alpha=alpha) +
        scale_x_continuous(breaks=c(.5,1,1.5), labels=c('0.5','1','1.5')) +
        # scale_y_continuous(breaks=c(0,.25,.5,.75,1), labels=c('0','','0.5','','1'), limits=c(0,1)) +
        labs(x=xlab, y=ylab, colour=collab) + 
        guides(colour=guide_legend(keyheight=.3, default.unit='inch'))
    p <- themefy(p)
}

Loading the data

allDataDir <- paste(db,'Projects/mc/data_bv3/mcBv3_xvv',sep='')
dataDirs <- dir(allDataDir)
dataDirs <- dataDirs[grep('xvv',dataDirs)]
colsOfInt <- c('participant', 'dom', 'session', 'mcBv', 'targTpeak', 'targXoff2', 
               'targV', 'stairStart', 'meanRev6')
df <- data.frame()
dfRevs <- data.frame()
dfIntn <- data.frame()
for(curDir in dataDirs){
    print(curDir)
    curDf <- read.csv(paste(allDataDir,'/',curDir,'/',curDir,'.csv',sep=''))
    curDf <- curDf[,colsOfInt]
    curDf$cond <- substr(curDir,19,22)
    subjStairs <- dir(paste(allDataDir,'/',curDir,'/',sep=''))
    subjStairs <- subjStairs[grep('.tsv',subjStairs)]
    for(curStairFN in subjStairs){
        curStair <- paste(allDataDir,'/',curDir,'/',curStairFN,sep='')
        curRevs <- read.table(curStair, skip=1, nrows=1)
        curIntn <- read.table(curStair, skip=4, nrows=2)
        curInfo <- readLines(curStair)
        curInfoCols <- data.frame(subj=curDf$participant[1], dom=curDf$dom[1], 
                                  sess=curDf$session[1], stairStart=curInfo[33], 
                                  mcBv=curInfo[41], targTpeak=curInfo[43],
                                  targEcc=curInfo[37], targV=curInfo[31])
        curDfRevs <- cbind(curInfoCols, curRevs[,2:9])
        nTrials <- ncol(curIntn)-1
        curDfIntn <- curInfoCols[rep(seq_len(nrow(curInfoCols)), each=nTrials),]
        curDfIntn$trial <- seq(1,nTrials)
        curDfIntn$intn <- as.numeric(curIntn[1,2:(nTrials+1)])
        curDfIntn$resp <- as.numeric(curIntn[2,2:(nTrials+1)])
        rownames(curDfIntn) <- NULL
        dfRevs <- rbind(dfRevs, curDfRevs)
        dfIntn <- rbind(dfIntn, curDfIntn)
    }
    df <- rbind(df, curDf)
}
## [1] "mc2_tgT-mcBv3_xvv-cent_p1_s1_2017-10-05_1028"
## [1] "mc2_tgT-mcBv3_xvv-cent_p12_s1_2017-10-10_1604"
## [1] "mc2_tgT-mcBv3_xvv-cent_p13_s1_2017-10-10_1653"
## [1] "mc2_tgT-mcBv3_xvv-cent_p5_s2_2017-10-06_1627"
## [1] "mc2_tgT-mcBv3_xvv-cent_p7_s2_2017-10-06_1432"
## [1] "mc2_tgT-mcBv3_xvv-dyna_p10_s1_2017-10-10_1033"
## [1] "mc2_tgT-mcBv3_xvv-dyna_p11_s2_2017-10-10_1425"
## [1] "mc2_tgT-mcBv3_xvv-dyna_p14_s2_2017-10-13_1137"
## [1] "mc2_tgT-mcBv3_xvv-dyna_p2_s1_2017-10-05_1234"
## [1] "mc2_tgT-mcBv3_xvv-dyna_p6_s1_2017-10-06_1040"
## [1] "mc2_tgT-mcBv3_xvv-dyna_p9_s1_2017-10-06_1549"
## [1] "mc2_tgT-mcBv3_xvv-peri_p13_s2_2017-10-12_1608"
## [1] "mc2_tgT-mcBv3_xvv-peri_p4_s1_2017-10-05_1539"
## [1] "mc2_tgT-mcBv3_xvv-peri_p5_s1_2017-10-05_1640"
## [1] "mc2_tgT-mcBv3_xvv-peri_p7_s1_2017-10-06_1245"
## [1] "mc2_tgT-mcBv3_xvv-stat_p11_s1_2017-10-10_1247"
## [1] "mc2_tgT-mcBv3_xvv-stat_p14_s1_2017-10-12_1240"
## [1] "mc2_tgT-mcBv3_xvv-stat_p2_s2_2017-10-05_1425"
## [1] "mc2_tgT-mcBv3_xvv-stat_p3_s1_2017-10-05_1342"
## [1] "mc2_tgT-mcBv3_xvv-stat_p6_s2_2017-10-12_1035"
## [1] "mc2_tgT-mcBv3_xvv-stat_p8_s1_2017-10-06_1357"
## [1] "mc2_tgT-mcBv3_xvv-stat_p9_s2_2017-10-10_1348"
ds <- rename(df, c(participant='subj', session='sess', meanRev6='thresh',
                   mcBv='maskV', targXoff2='targEcc'))
ds$targEcc <- round(ds$targEcc / 35,1)
ds$targV <- round(ds$targV / 3.5,1)
ds$maskType <- ''
ds$maskType[ds$maskV==0.01] <- 'static'
ds$maskType[ds$maskV==0.6] <- 'slow'
ds$maskType[ds$maskV==9.6] <- 'fast'
ds$maskV <- round(ds$maskV * 60 / 35, 1)
ds$maskV[ds$maskV<0.05] <- 0
dfIntn$maskV <- round(as.numeric(levels(dfIntn$mcBv))[dfIntn$mcBv] * 60 / 35, 1)
dfIntn$maskV[dfIntn$maskV<0.05] <- 0
ds$condSplit <- ''
ds$condSplit[ds$cond=='stat' | ds$cond=='dyna'] <- 'v'
ds$condSplit[ds$cond=='cent' | ds$cond=='peri'] <- 'ecc'
head(ds)
##   subj dom sess maskV targTpeak targEcc targV stairStart     thresh cond maskType condSplit
## 1    1   1    1   0.0       0.5     0.8   0.0         -2 -0.9583333 cent   static       ecc
## 2    1   1    1   1.0       0.5     0.8   4.6         -2 -0.7166667 cent     slow       ecc
## 3    1   1    1  16.5       1.0     0.8   0.0         -2 -1.6083333 cent     fast       ecc
## 4    1   1    1   1.0       1.5     0.8   4.6          0 -0.6583333 cent     slow       ecc
## 5    1   1    1   1.0       1.0     0.8   0.0         -2 -0.8750000 cent     slow       ecc
## 6    1   1    1  16.5       0.5     0.8   0.0         -2 -1.5916667 cent     fast       ecc

Transformed data sets

By target

thresh <- ddply(ds, .(subj,dom,sess,maskV,maskType,targTpeak,targEcc,targV,
                      cond,condSplit), summarise, threshMean = mean(thresh))
head(thresh)
##   subj dom sess maskV maskType targTpeak targEcc targV cond condSplit threshMean
## 1    1   1    1     0   static       0.5     0.8   0.0 cent       ecc -0.8916667
## 2    1   1    1     0   static       0.5     0.8   4.6 cent       ecc -1.2208333
## 3    1   1    1     0   static       1.0     0.8   0.0 cent       ecc -1.1041667
## 4    1   1    1     0   static       1.0     0.8   4.6 cent       ecc -1.2333333
## 5    1   1    1     0   static       1.5     0.8   0.0 cent       ecc -1.0875000
## 6    1   1    1     0   static       1.5     0.8   4.6 cent       ecc -1.2375000

Plotting across target eccentricities

Static target & static mask

ss <- sumFn(thresh[(thresh$maskType=='static' & thresh$targV==0),])
ssAve <- ss[ss$subj=='average',]
ssIndiv <- ss[ss$subj!='average',]
pTargStatMaskStat <- themefy(plotAve(ssAve)) + ylim(yLim)

Individual plots

(plotIndiv(ssIndiv) + ylim(-2,0))

Stair plots

iss <- dfIntn[dfIntn$maskV==0 & dfIntn$targV==0,]
(p <- ggplot(iss, aes(x=trial, y=intn, colour=targEcc, linetype=stairStart)) +
    facet_wrap(subj ~ targTpeak, ncol=6) + geom_point() + geom_line() + ylim(-2,0))

Dynamic target and static mask

ss <- sumFn(thresh[thresh$maskType=='static' & thresh$targV>0,])
ssAve <- ss[ss$subj=='average',]
ssIndiv <- ss[ss$subj!='average',]
pTargDynaMaskStat <- themefy(plotAve(ssAve)) + ylim(yLim)

Individual plots

(plotIndiv(ssIndiv) + ylim(-2,0))

Stair plots

iss <- dfIntn[dfIntn$maskV==0 & dfIntn$targV==16,]
(p <- ggplot(iss, aes(x=trial, y=intn, colour=targEcc, linetype=stairStart)) +
    facet_wrap(subj ~ targTpeak, ncol=6) + geom_point() + geom_line() + ylim(-2,0))

Static target & slow mask

ss <- sumFn(thresh[thresh$maskType=='slow' & thresh$targV==0,])
ssAve <- ss[ss$subj=='average',]
ssIndiv <- ss[ss$subj!='average',]
pTargStatMaskSlow <- themefy(plotAve(ssAve)) + ylim(yLim)

Individual plots

(plotIndiv(ssIndiv) + ylim(-2,0))

Stair plots

iss <- dfIntn[dfIntn$maskV==1 & dfIntn$targV==0,]
(p <- ggplot(iss, aes(x=trial, y=intn, colour=targEcc, linetype=stairStart)) +
    facet_wrap(subj ~ targTpeak, ncol=6) + geom_point() + geom_line() + ylim(-2,0))

Dynamic target and slow mask

ss <- sumFn(thresh[thresh$maskType=='slow' & thresh$targV>0,])
ssAve <- ss[ss$subj=='average',]
ssIndiv <- ss[ss$subj!='average',]
pTargDynaMaskSlow <- themefy(plotAve(ssAve)) + ylim(yLim)

Individual plots

(plotIndiv(ssIndiv) + ylim(-2,0))

Stair plots

iss <- dfIntn[dfIntn$maskV==1 & dfIntn$targV==16,]
(p <- ggplot(iss, aes(x=trial, y=intn, colour=targEcc, linetype=stairStart)) +
    facet_wrap(subj ~ targTpeak, ncol=6) + geom_point() + geom_line() + ylim(-2,0))

Static target & fast mask

ss <- sumFn(thresh[thresh$maskType=='fast' & thresh$targV==0,])
ssAve <- ss[ss$subj=='average',]
ssIndiv <- ss[ss$subj!='average',]
pTargStatMaskFast <- themefy(plotAve(ssAve)) + ylim(yLim)

Individual plots

(plotIndiv(ssIndiv) + ylim(-2,0))

Stair plots

iss <- dfIntn[dfIntn$maskV==16.5 & dfIntn$targV==0,]
(p <- ggplot(iss, aes(x=trial, y=intn, colour=targEcc, linetype=stairStart)) +
    facet_wrap(subj ~ targTpeak, ncol=6) + geom_point() + geom_line() + ylim(-2,0))

Dynamic target and fast mask

ss <- sumFn(thresh[thresh$maskType=='fast' & thresh$targV>0,])
ssAve <- ss[ss$subj=='average',]
ssIndiv <- ss[ss$subj!='average',]
pTargDynaMaskFast <- themefy(plotAve(ssAve)) + ylim(yLim)

Individual plots

(plotIndiv(ssIndiv) + ylim(-2,0))

Stair plots

iss <- dfIntn[dfIntn$maskV==16.5 & dfIntn$targV==16,]
(p <- ggplot(iss, aes(x=trial, y=intn, colour=targEcc, linetype=stairStart)) +
    facet_wrap(subj ~ targTpeak, ncol=6) + geom_point() + geom_line() + ylim(-2,0))

Group plots

grid.arrange(pTargStatMaskStat + #theme(legend.position='none') + 
                 ggtitle('a: static mask, static target'),
             pTargDynaMaskStat + ggtitle('d: static mask, dynamic target'),
             pTargStatMaskSlow + ggtitle('b: slow mask, static target'),
             pTargDynaMaskSlow + ggtitle('e: slow mask, dynamic target'),
             pTargStatMaskFast + ggtitle('c: fast mask, static target'),
             pTargDynaMaskFast + ggtitle('f: fast mask, dynamic target'),
             ncol=2, nrow=3) #widths=c((1-w)/2,w/2))

Analyses

Scaling function

cent <- function(v){
    v <- apply(v,2,function(x){
        x <- x - mean(unique(x),na.rm=T)
        x <- x / max(x)
    })
}

Centered data set

dsc <- ds
centCols <- c('dom','targTpeak','stairStart','targEcc','maskV','targV','sess')
dsc[,centCols] <- cent(dsc[,centCols])
# dsc$maskV <- (dsc$maskV / max(dsc$maskV)) #*2 - 1
dsc$maskV_c <- (dsc$maskV + 1) / 2
dsc$targV_c <- (dsc$targV + 1) / 2
dsc$targEcc_c <- round((dsc$targEcc + 1) / 2,0)
dsc$targTpeak_c <- (dsc$targTpeak + 1) / 2
# dsc$subj <- as.factor(dsc$subj)
# dsc$maskType <- ordered(ds$maskType, levels=c('static','slow','fast'))
dsc$maskType <- factor(dsc$maskType, c('static','slow','fast'))
head(dsc)
##   subj dom sess     maskV targTpeak targEcc targV stairStart     thresh cond maskType condSplit   maskV_c
## 1    1   1   -1 -0.546875        -1      -1    -1         -1 -0.9583333 cent   static       ecc 0.2265625
## 2    1   1   -1 -0.453125        -1      -1     1         -1 -0.7166667 cent     slow       ecc 0.2734375
## 3    1   1   -1  1.000000         0      -1    -1         -1 -1.6083333 cent     fast       ecc 1.0000000
## 4    1   1   -1 -0.453125         1      -1     1          1 -0.6583333 cent     slow       ecc 0.2734375
## 5    1   1   -1 -0.453125         0      -1    -1         -1 -0.8750000 cent     slow       ecc 0.2734375
## 6    1   1   -1  1.000000        -1      -1    -1         -1 -1.5916667 cent     fast       ecc 1.0000000
##   targV_c targEcc_c targTpeak_c
## 1       0         0         0.0
## 2       1         0         0.0
## 3       0         0         0.5
## 4       1         0         1.0
## 5       0         0         0.5
## 6       0         0         0.0

GLM

# ds$targTpeak_c <- ds$targTpeak - 0.5 # [.5, 1, 1.5] -> [0, .5, 1]
# ds$targEcc_c <- cent(ds$targEcc)
pvalfn(lmer(thresh ~ dom + stairStart + sess + maskType * targTpeak_c * targEcc_c +
                 (1|subj), data=dsc[dsc$targV==-1,]))
##                                      estm    se      df    tVal  pVal star
## (Intercept)                        -1.027 0.063  22.260 -16.315 0.000  ***
## dom                                 0.005 0.053  10.828   0.092 0.928     
## stairStart                          0.030 0.009 388.056   3.150 0.002   **
## sess                                0.006 0.018 355.291   0.340 0.734     
## maskTypeslow                        0.184 0.051 388.056   3.646 0.000  ***
## maskTypefast                       -0.325 0.051 388.056  -6.421 0.000  ***
## targTpeak_c                        -0.071 0.055 388.056  -1.278 0.202     
## targEcc_c                           0.068 0.052 389.444   1.308 0.192     
## maskTypeslow:targTpeak_c            0.027 0.078 388.056   0.345 0.730     
## maskTypefast:targTpeak_c            0.016 0.078 388.056   0.199 0.842     
## maskTypeslow:targEcc_c             -0.310 0.073 388.056  -4.231 0.000  ***
## maskTypefast:targEcc_c             -0.093 0.073 388.056  -1.271 0.205     
## targTpeak_c:targEcc_c               0.010 0.080 388.056   0.128 0.899     
## maskTypeslow:targTpeak_c:targEcc_c  0.042 0.113 388.056   0.369 0.712     
## maskTypefast:targTpeak_c:targEcc_c  0.025 0.113 388.056   0.216 0.829
pvalfn(lmer(thresh ~ dom + stairStart + sess + maskType * targTpeak_c * targEcc_c *
                 targV_c + (1|subj), data=dsc))
##                                              estm    se      df    tVal  pVal star
## (Intercept)                                -1.075 0.068  20.898 -15.905 0.000  ***
## dom                                         0.036 0.059  11.908   0.611 0.552     
## stairStart                                  0.029 0.006 752.920   4.465 0.000  ***
## sess                                       -0.027 0.008 758.904  -3.538 0.000  ***
## maskTypeslow                                0.184 0.048 752.920   3.839 0.000  ***
## maskTypefast                               -0.325 0.048 752.920  -6.761 0.000  ***
## targTpeak_c                                -0.071 0.053 752.920  -1.345 0.179     
## targEcc_c                                   0.062 0.049 753.271   1.262 0.207     
## targV_c                                     0.054 0.049 753.263   1.095 0.274     
## maskTypeslow:targTpeak_c                    0.027 0.074 752.920   0.364 0.716     
## maskTypefast:targTpeak_c                    0.016 0.074 752.920   0.210 0.834     
## maskTypeslow:targEcc_c                     -0.310 0.069 752.920  -4.455 0.000  ***
## maskTypefast:targEcc_c                     -0.093 0.069 752.920  -1.338 0.181     
## targTpeak_c:targEcc_c                       0.010 0.076 752.920   0.134 0.893     
## maskTypeslow:targV_c                        0.054 0.069 752.920   0.770 0.442     
## maskTypefast:targV_c                        0.465 0.069 752.920   6.696 0.000  ***
## targTpeak_c:targV_c                         0.035 0.076 752.920   0.458 0.647     
## targEcc_c:targV_c                           0.048 0.071 752.923   0.677 0.499     
## maskTypeslow:targTpeak_c:targEcc_c          0.042 0.108 752.920   0.389 0.698     
## maskTypefast:targTpeak_c:targEcc_c          0.025 0.108 752.920   0.228 0.820     
## maskTypeslow:targTpeak_c:targV_c           -0.018 0.108 752.920  -0.167 0.867     
## maskTypefast:targTpeak_c:targV_c            0.056 0.108 752.920   0.516 0.606     
## maskTypeslow:targEcc_c:targV_c              0.120 0.101 752.920   1.196 0.232     
## maskTypefast:targEcc_c:targV_c             -0.009 0.101 752.920  -0.086 0.931     
## targTpeak_c:targEcc_c:targV_c              -0.065 0.110 752.920  -0.587 0.558     
## maskTypeslow:targTpeak_c:targEcc_c:targV_c  0.072 0.156 752.920   0.462 0.644     
## maskTypefast:targTpeak_c:targEcc_c:targV_c -0.007 0.156 752.920  -0.042 0.966